clear
set virtual off
set more off
pause off
clear all
 
local vecY 10 17

local doestimation=1
local dosigma=0

if `doestimation'==1 {

use processed/UKquasi_panel.dta , clear

 
local inityear=1975
local finalyear=2017

local pickdiff=2

keep if year>=`inityear' & year<=`finalyear'

 foreach n of numlist 1/17 {
 
	gen logprice`n'=log(price`n')  
 }


reshape long logprice b_, i(year pct1000) j(good)

rename b_ bshare

gen chooseY=0
foreach x of numlist `vecY' { 
	replace chooseY=1 if good==`x'
}

gen chooseX=1-chooseY

sort pct1000 year good

bys pct1000 year: egen bX=sum(bshare) if chooseX==1

keep if chooseX==1

keep year pct1000 good totexp_mean totexp_median logprice bshare bX

gen bXi=bshare/bX
gen logbX=log(bX)


//merge m:1 year good using Engel_coef_EDIT.dta

drop if pct1000==.

gen logoil=.

replace logoil=1.489 if year==1974
replace logoil=1.614 if year==1975
replace logoil=1.947 if year==1976
replace logoil=2.103 if year==1977
replace logoil=2.046 if year==1978
replace logoil=2.357 if year==1979
replace logoil=2.777 if year==1980
replace logoil=2.897 if year==1981
replace logoil=2.957 if year==1982
replace logoil=2.998 if year==1983
replace logoil=3.086 if year==1984
replace logoil=3.071 if year==1985
replace logoil=2.327 if year==1986
replace logoil=2.458 if year==1987
replace logoil=2.193 if year==1988
replace logoil=2.481 if year==1989
replace logoil=2.619 if year==1990
replace logoil=2.498 if year==1991
replace logoil=2.455 if year==1992
replace logoil=2.509 if year==1993
replace logoil=2.418 if year==1994
replace logoil=2.457 if year==1995
replace logoil=2.653 if year==1996
replace logoil=2.532 if year==1997
replace logoil=2.161 if year==1998
replace logoil=2.477 if year==1999
replace logoil=2.995 if year==2000
replace logoil=2.891 if year==2001
replace logoil=2.855 if year==2002
replace logoil=2.947 if year==2003
replace logoil=3.118 if year==2004
replace logoil=3.435 if year==2005
replace logoil=3.580 if year==2006
replace logoil=3.588 if year==2007
replace logoil=3.983 if year==2008
replace logoil=3.674 if year==2009
replace logoil=3.940 if year==2010
replace logoil=4.082 if year==2011
replace logoil=4.085 if year==2012
replace logoil=4.137 if year==2013
replace logoil=4.036 if year==2014
replace logoil=3.461 if year==2015
replace logoil=3.460 if year==2016
replace logoil=3.676 if year==2017

sort pct1000 good year

bys pct1000 good: gen dlogpoil=logoil-logoil[_n-`pickdiff'] if year-year[_n-`pickdiff']==`pickdiff'


//reghdfe logprice [aweight=bshare], absorb(good year) res(residlogp)
//reghdfe logprice, absorb(good year) res(residlogp)


//bys pct1000 good: gen dresidlogp=residlogp-residlogp[_n-1] if year-year[_n-1]==1
bys pct1000 good: gen dresidlogp=logprice-logprice[_n-`pickdiff'] if year-year[_n-`pickdiff']==`pickdiff'

bys pct1000 good: gen dlogbX=logbX-logbX[_n-`pickdiff'] if year-year[_n-`pickdiff']==`pickdiff'

gen sh_dlogp=dresidlogp*bXi

bys pct1000 year: egen sumshares=sum(bXi)
bys pct1000 year: egen dlogp_pct=sum(sh_dlogp)

//keep if year>`inityear'


collapse dlogp_pct dlogbX dlogpoil bX, by (year pct1000) 

gen indH=0
replace indH=1 if pct1000>500
gen dlogpH=indH*dlogp_pct
gen dlogpoilH=indH*dlogpoil



}




* Regressions
reghdfe dlogbX dlogp_pct, cluster(pct1000) absorb(year pct1000)
estimates store OLS1
ivreghdfe dlogbX (dlogp_pct=dlogpoil), cluster(pct1000) absorb(pct1000)
estimates store TSLS

ivreghdfe dlogbX (dlogp_pct=dlogpoil) , cluster(pct1000) absorb(pct1000) first savefirst savefprefix(st1)
estadd scalar FF = `e(widstat)' : TSLS

reghdfe dlogbX dlogp_pct dlogpH , cluster(pct1000) absorb(year pct1000)
estimates store OLS2
ivreghdfe dlogbX (dlogp_pct dlogpH =dlogpoil dlogpoilH), cluster(pct1000) absorb(pct1000)
estimates store TSLS2
ivreghdfe dlogbX (dlogp_pct dlogpH =dlogpoil dlogpoilH) , cluster(pct1000) absorb(pct1000) first savefirst savefprefix(st1)
estadd scalar FF = `e(widstat)' : TSLS2


estfe . OLS1 TSLS OLS2 TSLS2, labels(year "Year" pct1000 "HH")

esttab OLS1 TSLS OLS2 TSLS2 , se b(3) star(* 0.1 ** 0.05 *** 0.01) nodepvar keep(dlogp_pct dlogpH) indicate(`r(indicate_fe)' ) ///
 coeflabels(dlogp_pct "dlogp" dlogpH "dlogp $\times$ $\mathbb{I}_{nt}^{H}$") label replace stats(FF FE N , fmt(0 0)  labels(`"F-stat"' `"FE"' `"Obs"')) scalar(FF) ///
   nonotes ///
collabels("Model 1" "Model 2" "Model 3" "Model 4" "Model 5", none) 



